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We show that multichannel quantum defect theory (MQDT) can be applied successfully as an 
efficient computational method for cold molecular collisions in Li+NH, which has a deep and strongly 
anisotropic interaction potential. In this strongly coupled system, closed-channel poles restrict the 
range over which the MQDT Y can be interpolated. We present an improved procedure to transform 
the MQDT reference functions so that the poles are removed from the energy range of interest. 
Effects due to very long-range spin-dipolar couplings are outside the scope of MQDT, but can be 
added perturbatively. The new procedure makes it possible to calculate the elastic and inelastic cross 
sections needed to evaluate the feasibility of sympathetic cooling of NH by Li using coupled-channel 
calculations at only 5 combinations of energy and field. 



I. INTRODUCTION 

Cold and ultracold molecules provide an exciting door- 
way to new fields in physics and chemistry. They share 
the controllability and tunable interactions that have 
made cold atom studies so fruitful. However, the richer 
structure of molecules makes them suitable for many 
new applications and research directions [TJ [2]. Polar 
molecules are of particular interest because their electric 
dipole moment allows the interactions between them to 
be controlled and manipulated by external fields. This 
strong tunable response may make it possible to develop 
a fully controlled chemistry [3 where every degree of free- 
dom of the reaction can be tuned, providing fundamental 
insights into chemical reaction processes. Molecules such 
as KRb [4] and Cs2 [5] have already been produced at 
submicrokelvin temperatures in their lowest-energy elec- 
tronic, vibrational, rotational, and hyperfme state, by 
magnetoassociation followed by laser state transfer. Re- 
actions of ultracold 40 K 87 Rb with itself and with K and 
Rb atoms have been studied [6], and it was seen that 
quantum statistics and quantum threshold laws play an 
important role in determining the rates of inelastic colli- 
sions. 

The only ultracold molecules that are accessible with 
current methods are alkali-metal dimers. However, there 
is great interest in producing samples of other molecu- 
lar species in the ultracold regime. The most promising 
route to this is first to cool and trap them in the cold 
regime (at temperatures below 1 K), using a method such 
as buffer-gas cooling [7] or molecular beam deceleration 
[8] and then to bring them to the ultracold regime using 
a second-stage approach such as evaporative cooling [9], 
sympathetic cooling [10] , or laser cooling [TTJ [12] . 

Evaporative and sympathetic cooling both rely on elas- 
tic collisions to thermalize the sample, but both can be 
prevented by inelastic collisions that release kinetic en- 
ergy and lead to trap loss. Collisional properties are 
also key to methods for controlling ultracold atoms and 
molecules. Calculations on atomic and molecular colli- 



sions are therefore crucial to both the production and 
control of cold and ultracold molecules. Such calcula- 
tions require the solution of the set of coupled differ- 
ential equations obtained from the Schrodinger equation. 
There are are various numerical methods for solving these 
coupled differential equations, of which the most com- 
monly used is the full coupled-channels method. This 
propagates the matrix solution of the Schrodinger equa- 
tion from short to long range and takes a time propor- 
tional to A 7 " 3 , where N is the number of coupled channels. 
The properties of completed collisions are described by 
the scattering matrix S, which is obtained by matching 
the propagated solutions to free-particle wavefunctions 
(Ricatti-Bessel functions) at long range [T3] . 

Full coupled-channels calculations can be extremely 
expensive, particularly in applied electric and magnetic 
fields [T4HT6] . The expense is particularly great for sys- 
tems with deep and strongly anisotropic potential wells 
[T7] , for molecule- molecule collisions [18] , or when nuclear 
hyperfme interactions are included p~9j [20] . 

In cold collision studies, the scattering S matrix is of- 
ten a fast function of collision energy E and magnetic 
field £>, with extensive structure due to scattering res- 
onances and discontinuous behavior at threshold. Cal- 
culations are thus required over a fine grid of energies 
and/or applied field, typically over an energy range of 
order 1 K, from threshold, and for magnetic fields up to 
a few thousand G. This contrasts with the situation for 
collisions of ultracold atoms, where the energy range of 
interest is commonly a few fiK and the fields are typically 
a few hundred G. 

Approaches based on multichannel quantum defect 
theory (MQDT) avoid the repetition of the expensive 
propagation by representing the scattering properties in 
terms of a matrix Y(E, B) that is a smooth function of E 
and B [2TH25], MQDT has proved highly successful for 
cold atomic interactions [24-30 , and more recently it has 
been applied to collisions of cold and ultracold molecules 
[3TH36] . MQDT defines the matrix Y(E, B) at a match- 
ing distance, r matc h, at relatively short range in order to 
achieve this insensitivity to energy and field. The Y ma- 
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trix contains all the scattering dynamics inside r matc h- 
The smooth variation of Y allows it to be obtained once 
and then used for calculations over a wide range of ener- 
gies and fields, or obtained by interpolation from a few 
points. The computational cost of calculations at addi- 
tional energies and fields is only proportional to N, not 
N 3 . 

We have previously demonstrated the application of 
MQDT to cold molecular collisions for the moderately 
anisotropic system Mg+NH( 3 S-) PI [36]. In ref. [55]. we 
showed that the choice made for the phase of the MQDT 
reference functions is very important in producing a Y 
matrix that can be interpolated smoothly over a wide 
range of energy and field. The purpose of the present pa- 
per is to explore how the approach performs for a much 
more strongly anisotropic system with many more closed 
channels. For this purpose we choose Li+NH, which has 
been studied previously using full coupled-channels cal- 
culations by Wallis et al. [17] . As in ref. [17] , we focus on 
collisions between spin-polarized Li and NH, which oc- 
cur on the quartet potential energy surface. This surface 
is deep and highly anisotropic, with a well depth about 
1800 cm -1 at the Li-NH geometry but only 113 cm -1 at 
the NH-Li geometry. With a small but important mod- 
ification to the method of choosing phases described in 
ref. [36] , we find that we can obtain accurate results for 
elastic and inelastic cross section, over the entire range 
of energies relevant to sympathetic cooling, using only 5 
coupled-channels calculations. 

II. THEORY 



over a wide range of energy and field [36]. Equation ([!]) 
shows that the Y matrix has a pole whenever the com- 
ponent of the propagated multichannel wavefunction ipi 
in any channel i is proportional to the reference func- 
tion gi and has no contribution from fa, i.e. when gi and 
the full coupled-channels solution have the same phase at 
r matc h- In r ef. [36 , we proposed transforming the pair of 
reference functions fa and gi in channel i with a rotation 
angle 0$, chosen so that the diagonal matrix elements Ya 
is at a particular reference energy E re f and field B re f. 
This ensures that the reference function ^ and the full 
coupled-channels solution in channel i are perfectly out 
of phase at the chosen r matc h, and the resulting Y matrix 
is therefore pole- free close to E re f and B ve t [36] . 

The range of the pole-free region is dependent on where 
the matching occurs. When matching is in the classically 
allowed region, the phases of both the reference functions 
and the propagated coupled-channels solutions vary ap- 
proximately linearly with energy and setting the diagonal 
Y matrix elements to zero is effective: the relative phase 
of the reference functions and the coupled-channels so- 
lution is a slow function of energy. For a closed channel 
where matching is carried out in the classically forbidden 
region, however, there is resonance structure in both the 
coupled-channels solutions and the reference functions so 
that the phase of each is a fast (and nonlinear) function 
of energy. Fortunately, the energies at which the reso- 
nance structure occurs depend on the choice of d{. In the 
present work we show that a more sophisticated choice 
of 6i than that of ref. [36 can produce a larger pole-free 
region for closed channels. 



The theory of MQDT is given in details in refs. [22H25] . 
Here we give only a brief description, following references 
l34l l36l , which is sufficient to describe the notation we 



use. 



MQDT makes the approximation that the multichan- 
nel Schrodinger equation is uncoupled at distances r > 
Rmatch, so that its solution in this region may be written 
in the matrix form 



(1) 



Here / and g are diagonal matrices containing the func- 
tions fa and gi, which are linearly independent solutions 
of a reference Schrodinger equation in each asymptotic 
channel z, 



' 2[i dr 2 



+ Uf{r) 



Mr) = 0, 



(2) 



and similarly for gi(r). ji is the reduced mass, and 
the reference potentials E/J ef (r) approach the true po- 
tential at long range. They include the centrifugal terms 
h 2 Li(Li + l)/2/ir 2 , where Li is the partial-wave quantum 
number for channel i. 

The phase of the short-range reference functions fa and 
gi is a disposable parameter of MQDT and may be cho- 
sen to generate a Y matrix that is smooth and pole-free 



A. Basis set and quantum numbers 

We construct the collision Hamiltonian in the fully 
uncoupled basis set |nra n )|sNH^s N H)l 5 Li^s Li )|£^L), 
where the quantum numbers n and snh describe the 
rotation and electron spin of the NH molecule and 
describes the electron spin of the Li atom. The corre- 
sponding ra quantum numbers are the projections onto 
the space-fixed magnetic field axis. Hyperfine structure 
is neglected. The matrix elements required for the cou- 
pled equations are the same as for scattering of NH 
from a closed-shell atom [16], with the addition of the 
anisotropic intermolecular spin-spin dipolar interaction 

mi. 

The coupled-channels calculations may in principle be 
carried out in any sufficiently complete basis set. How- 
ever, the Y and S matrices are expressed in a basis set 
of eigenfunctions of the field-dressed monomer Hamil- 
tonians. At low field the states of the NH molecule 
are approximately described by quantum numbers j and 
raj, where j is the resultant of n and s. We label ele- 
ments of Y and S by subscripts n, j, ra^, ra SLi , L, Ml — » 
??/, j', ra^, ra^ L ., L', M' L . For diagonal elements we sup- 
press the second set of labels. 



3 



Numerical methods 



The coupled-channels calculations required for both 
MQDT and the full coupled-channels approach are car- 
ried out using the MOLSCAT package [37], as modified 
to handle collisions in magnetic fields [16]. The cou- 
pled equations are solved numerically using the hybrid 
log-derivative propagator of Alexander and Manolopou- 
los [38], which uses a fixed-step-size log-derivative prop- 
agator in the short-range region (r m i n < r < r m id) and a 
variable-step-size Airy propagator in the long-range re- 
gion (r m id < r < r max ). As in ref. [17 , the full coupled- 
channels calculations use r m i n = 1.8 A, r mid = 12.5 A 
and r max = 600 A (where 1 A = lO" 10 m). MQDT 
requires coupled-channels calculations only from r m i n to 
r matc h (which is less than r m id), so only the fixed- step-size 
propagator is used in this case. 

The initial MQDT reference functions and quantum 
defect parameters are obtained as described in ref. [34] . 
using the renormalized Numerov method [39 to solve 
the 1-dimension Schrodinger equations for the reference 
potentials. In this paper all MQDT calculations use the 
reference potential 



Uf(r) = V (r) + 



h 2 U{U + 1) 
2/ir 2 



E? 



(3) 



where Vo(r) is the isotropic part of the interaction po- 
tential and Li is the partial-wave quantum number for 
channel i. The reference potential contains a hard wall 



at r 



wall 



so that U- ef (r) = oo for r < r™ a11 . In the 
present paper we choose r™ a11 = 4.0 A. This reference 
potential has been shown to produce quantitatively ac- 
curate results for Mg+NH [341136] . 

The uncoupled basis functions used to solve the 
coupled-channel equations are not eigenfunctions of the 
Hamiltonians of the separated monomers. The log- 
derivative matrix obtained from the coupled-channel cal- 
culations at a distance r matc h is therefore transformed 
into a basis set that diagonalizes the asymptotic Hamil- 
tonian [40]. The MQDT Y matrix is then obtained by 
matching to this log-derivative matrix at r matc h using Eq. 
(13) of ref. [34]. All channels with n > 2 are treated as 
strongly closed and thus not included in the MQDT part 
of the calculation, but are included in the log- derivative 
propagation. 



III. RESULTS AND DISCUSSION 

In a magnetic field, the lowest Li-NH threshold (n = 
0, snh = l^Li = ^) splits into 6 Zeeman sublevels. We 
consider collisions between Li atoms and NH molecules 
that are both initially in their magnetically trappable 
low-spin-seeking states, m Sy u = +1/2 and m S) NH = +1- 
This corresponds to the highest of the 6 thresholds. 

Figure [T] shows the variation of the representative el- 
ement ^i 5 2,-2,-^,7,45 as a function of the matching dis- 
tance and energy, when the phases 0i are chosen to make 




FIG. 1. (Color online) arctanl 7 ^ / V for a single representative 
diagonal Y matrix element, as a function of collision energy 
E and r ma tch, with the phases 0i set so that Yu = in all 
channels for energy E re { = 0.5 K and field B re { = 10 G. 



all diagonal Y matrix elements zero at collision energy 
Eref — 0.5 K and field B re f — 10 G. It may be seen that 
there are poles (where arctan Y^/V passes through ±1/2) 
whose positions depend strongly on r matc h- Other Y ma- 
trix elements are quantitatively different but have poles 
in the same places. The basis set used for Figure [I] in- 
cludes all functions up to n max = 3 and L max = 8. This 
unconverged basis set was used due to the substantial 
computational cost of performing a full coupled-channels 
calculation at every energy in the Figure. The outer turn- 
ing point of the n = 1 reference potential is at 6.1 A, and 
it may be seen that, for values of r matc h inside this, the Y 
matrix is pole- free over many K. However, MQDT with 
such small values of r matc h does not produce accurate re- 
sults because it neglects all channel couplings that exist 
outside r matc h- When r matc h > 6.1 A, poles start to enter 
the Y matrix in the energy range of interest. As r matc h 
increases further, the poles move and at some values of 
Rmatch can come within 0.1 K of E re f. 

A contour plot such as Fig. [I] requires coupled-channels 
calculations at every energy, and producing it thus sac- 
rifices most of the computational savings that MQDT is 
designed to achieve. In addition, we need a procedure for 
choosing the phases 6i that will guarantee a large pole- 
free region for any choice of r matc h- In the remainder 
of this paper, we perform calculations at only a single 
value of r matc h = 6.5 A, deliberately chosen to be in a 
region where Fig. [I] shows that there are poles in Y in- 
conveniently close to E Ye {. In addition, the remaining 
calculations use a basis set including all functions up to 
n max = 6 and L max = 8, except where stated otherwise. 

The dependence of a diagonal Y matrix element on the 
phase of the reference functions (in any channel handled 
by MQDT) is 



Y u {E)=ttm [0i + 5i(E)], 



(4) 
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FIG. 2. (Color online) Upper panel: arctan / 'ir from 
coupled-channels calculations at B = 10 G for a single di- 
agonal element, ^2, -2, -±,7,45 w ^ n ®%' se t so that Y^v = in 
all other channels z' at E Te { = 0.5 K (horizonal dashed line). 
Lower panel: arctan Y^/7r as a function of energy and 6i for 
a single uncoupled channel as given by Eqs. Q and |5]). 

where Si(E) is the phase shift between the unrotated 
reference function fa and the propagated multichannel 
wavefunction in channel i. For a closed channel that is 
capable of supporting resonances, the phase shift around 
a resonance has a Breit-Wigner form, 

Si(E) = Si(E) + arctan (-3^) , (5) 

where Si is a slowly varying (non-resonant) background 
term. The resonant part of this function is shown in 
the lower panel of Fig. |2j for values of the parameters 
appropriate to one of the channels in Li+NH. It may 
be seen that choosing a value of 0i that makes Ya zero 
(shown by the dashed vertical line) does not guarantee 
a large pole-free region in the case where E re f is close to 
E\ es . A much better choice in this case is to set 0^ to the 
value shown by the solid vertical line. In the following 



we will show how this can be achieved. 

A basic problem of MQDT in coupled-channel prob- 
lems is that a pole in Y that originates in any channel 
causes a pole in every channel. We refer to this as the 
contamination of one channel by another. The upper 
panel of Figure [2] shows arctan Yh/tt for the single matrix 
element, ^1,2,-2,-^,7,45 obtained from coupled-channels 
calculations as a function of 0^ and collision energy E. 
The phases Oy in all other channels %' are set to the val- 
ues that produce Yyy = at E re ? = 0.5 K. The broad 
horizontal sweep around K arises from a resonance in 
channel z, while the narrower sweeps at —1.8 K, —1.6 K, 
— 1.2 K 1.4 K and 2.0 K are poles due to contamination 
from other channels. Setting Yyy to zero in all these 
other channels has shifted these contamination effects to 
energies either above about 1.4 K or below —1.2 K , leav- 
ing a region of about 2.6 K uncontaminated by other 
channels. For the specific circumstances shown in Figure 
|2j it is seen that choosing 0^ so that Ya = results in 
Oi/n « 0.33 (dashed vertical line). This produces a pole 
in Ya itself at about —0.1 K, which is inconveniently close 
to E xe {. However, a choice of Oi/n = 0.63 (solid vertical 
line) would produce a much larger pole-free range, lim- 
ited only by poles due to contamination effects in other 
channels. If improved values of 0y can also be obtained 
for these contamination poles, then there is clearly the 
prospect of achieving a much improved pole-free region. 

The pole structure in channel i when uncontaminated 
by other channels is given by Eqs. Q and ([5]). In order to 
use these to obtain a better choice of 0^, we need values 
for E\ es , and Si. To obtain these we first optimize the 
phases as in ref . [36] , transforming the reference functions 
so that Ya = in all channels at energy E re f. This pro- 
vides at least a small region where Ya is uncontaminated 
by poles in other channels. We then carry out coupled- 
channels calculations at 2 additional energies near -E re f, 
and use Eqs. Q and ([5| to obtain the three parameters 
i^ es , and Si numerically, neglecting the slow varia- 
tion of Si with E. The optimum pole- free region for this 
channel is then achieved by setting 6i = tt/2 — Si. 

The pole-free region for the entire Y matrix is opti- 
mized by applying this procedure in all channels where 
there is resonant structure close to the reference energy. 
We first calculate the numerical second derivative of the 
diagonal matrix elements Ya with respect to energy. We 
then select the channel with the largest second derivative, 
apply the procedure described above, and use the new set 
of phases to recalculate the three Y matrices. Because of 
channel mixing, this in principle changes all the diagonal 
matrix elements. If it reduces \d 2 Ya/dE 2 \ then we ac- 
cept the new value of 0^. If not, we move on to the next 
channel and apply the same procedure. We loop over 
the channels in this manner until there is no channel for 
which changing 0i to tt/2 — Si reduces \d 2 Yu/dE 2 \. 
This is an inexpensive procedure, as it uses the same 3 
coupled-channels calculations as before. Only the closed 
channels need to be included in the loop since only these 
channels have resonance structure. 



FIG. 3. (Color online) Yu as a function of energy for all 
channels, with all phases Oi optimized as described in ref. [3l] 
(dashed red lines) or using Eqs. Q and ^ (solid black lines). 
In both cases the phases are optimized at r ma tch = 6.5 A for 
E Tei = 0.5 K and B rei = 10 G. 

Figure [3] compares the final matrix elements Ya in all 
the channels included in the MQDT procedure, obtained 
with the two optimization schemes. The dashed red lines 
show the result of choosing Oi so that Ya is zero in every 
channel, as in ref. [36 , while the solid black lines show the 
result of optimizing Oi as described above. Both calcula- 
tions use r matc h = 6.5 A and optimize Oi at E re f — 0.5 K 
and B ve { = 10 G. It may be seen that taking account 
of closed-channel resonances significantly increases the 
pole- free range of Y. Furthermore, it produces Y ma- 
trix elements that are considerably more linear between 
and 1 K and may thus be interpolated more accurately. 

Figure [Z] compares diagonal T-matrix elements \Ta\ 2 
(where = 5^ — Sij) obtained from full coupled- 
channels calculations with those from the MQDT method 
using interpolation. The MQDT results were obtained by 
interpolating (and extrapolating) Y quadratically using 
3 points separated by 0.1 K around 0.5 K. The MQDT 
results obtained by interpolation are very similar to the 
full coupled-channels results, even around the resonance 
features at E w 0.7 K. 



IV. APPLICATION TO SYMPATHETIC 
COOLING 

The key quantity that determines whether sympathetic 
cooling can succeed is the ratio 7 of elastic to inelas- 
tic cross sections for collisions of trapped molecules with 
coolant atoms. This ratio typically needs to be greater 
than about 100 if trapped molecules are to undergo 
enough elastic (thermalizing) collisions to achieve cool- 
ing before undergoing an inelastic collision that releases 
kinetic energy and causes trap loss. Wallis et al. [17] have 
investigated Li + NH collisions using coupled-channel 



FIG. 4. (Color online) Squares of diagonal T-matrix elements 
T mj ,L,ML m the incoming channels for rrij = +1 and L = 0, 
2 and 6 at B = 10 G, obtained from full coupled-channels 
calculations (solid, black) and MQDT using optimized refer- 
ence functions for r ma tch = 6.5 A both with (dot-dash, blue) 
and without (dashed, red) interpolation. L = 4 is not shown 
because it obscures the resonant feature for L = 6. 



calculations, and produced contour plots that show the 
ratio 7 as a function of collision energy E and magnetic 
field B. These calculations were very expensive because 
they required calculations on a fine grid of energies and 
fields. The contour plots given in ref. [17] actually in- 
cluded coupled-channel calculations at 204 combinations 
of energy and field. 

MQDT offers the possibility of producing the entire 
contour plot from coupled-channels calculations at just 
three energies (required to optimize the phases) and two 
magnetic fields. We have therefore used MQDT to repeat 
the calculations of ref. [17 on the unsealed potential en- 
ergy surface. Coupled-channel calculations from r m i n to 
Rmatch were carried out at collision energies E = 0.01, 
0.05 and 0.1 K t magnetic field B = 10 G to optimize 
the phases, and at E = 0.01 and 0.1 K at B = 1000 G to 
allow linear interpolation in B and E. These calculations 
used a basis set with n max = 6 and L max = 8 to allow 
direct comparison with ref. [17]. The resulting contour 
plot of 7 is compared with the results of ref. [17] in the 
upper panel of Figure [5] it may be seen that there is good 
agreement at both high collision energies (E > 0.01 K) 
and high fields (B > 100 G), but that MQDT by itself 
breaks down when both E and B are small [4T] . 

The inaccuracy in MQDT at low energy and low field 
occurs because, in this region, the inelastic cross sections 
are dominated by long-range inelasticity involving the 
magnetic dipole interaction between the spins of Li and 
NH. As described by Janssen et al. [18] . there are long- 
range avoided crossings between the effective potential 
curves for the incoming channel and for inelastic chan- 
nels with larger values of L. These crossings usually occur 
outside the centrifugal barriers, and even coupled-channel 
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FIG. 5. (Color online) Contour plots of the ratio 7 of elas- 
tic to inelastic cross sections for Li+NH collisions as a func- 
tion of energy and magnetic field, obtained using a basis set 
with n max = 6 and L max = 8. The solid contours and shad- 
ing show the MQDT results and the dashed contours show 
the results from full coupled-channel calculations [T7]. The 
diagonal blue lines show the field at which the Zeeman en- 
ergy in a quadrupole trap is Upper panel: results 
from MQDT alone. Lower panel: results from MQDT with 
long-range spin-spin couplings added using the hybrid Born 
approach. 



interactions. The DWBA, by contrast, is quite accurate 
for initial L = but can be unstable when L and L' are 
both non-zero. However, in the latter case, short-range 
effects are unimportant. We have therefore used a hy- 
brid Born approximation, made up of the DWBA for ini- 
tial L = and the simple Born approximation for initial 
L > 0. When we add the resulting inelastic cross sections 
to those from MQDT, we obtain the contour plot for the 
ratio 7 shown in the lower panel of Figure [5] It may be 
seen that this gives essentially complete agreement with 
the full coupled-channel results. 

The MQDT approach makes it feasible to use a larger 
basis than was possible in ref. [17 and to carry out the 
calculations on a much finer grid of energies and fields. 
Figure [6] shows the results obtained from MQDT with 
perturbative long-range corrections for a converged basis 
set with n max = 10 and L max = 8, with cross sections 
calculated on a 51 x 51 grid of energies and fields. The 
resulting coupled-channels basis sets contains 1887 basis 
functions, as compared to 937 functions for the smaller 
basis set used in Figure |5j and each coupled-channel cal- 
culation is therefore a factor of 8 more expensive. 

As described in ref. [17] , the elastic and inelastic cross 
sections are a strong function of both potential scaling 
and basis set size, because they depend sensitively on 
the positions of near-dissociation levels. Because of this, 
calculations on a single potential do not give quanti- 
tative predictions for the ratio of elastic and inelastic 
cross sections, and it is essential to explore the potential- 
dependence of the results. We found that using the un- 
sealed potential with a converged basis set gave a highly 
atypical contour plot, because it has an accidentally near- 
zero scattering length and therefore a very small elastic 
cross section. The calculations in Figure [6] used a poten- 
tial with an overall scaling factor of 0.995, which produces 
a much more typical contour plot. Exploring the depen- 
dence of the results on the scaling factor confirmed the 
conclusions of ref. [17], that sympathetic cooling of NH 
by Li is likely to succeed if the molecules can be precooled 
to a temperature around 20 mK. 



V. CONCLUSIONS 



calculations must be propagated to very long range (hun- 
dreds of A) to capture their effects. They are thus outside 
the scope of MQDT, which neglects couplings outside 

^match • 

The long-range couplings may however be included 
perturbatively at very little expense. Janssen et al [42] 
have developed both a simple Born approximation and a 
distorted-wave Born approximation (DWBA) for calcu- 
lating the inelastic cross sections due to these long-range 
couplings. The simple Born approximation is stable to 
evaluate, but can give results up to a factor of 2 in er- 
ror for cross sections for initial L = because it does 
not take account of the phase shift due to short-range 



We have demonstrated that Multichannel Quantum 
Defect Theory (MQDT) can provide quantitatively ac- 
curate cross sections for cold and ultracold elastic and 
inelastic collisions in magnetic fields for a strongly cou- 
pled molecular system, Li + NH. However, the choice of 
the phases of the MQDT reference functions is crucial. 
For Mg+NH, it was sufficient to choose the phases so that 
Ya = in every channel included in MQDT [36 . For 
Li+NH, however, this does not guarantee that all closed- 
channel poles are moved far away in energy, and the poles 
can cause problems in interpolation. In the present pa- 
per, we have developed an improved approach for opti- 
mizing the phases that ensures that closed-channel poles 
are far away from a reference energy. 
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FIG. 6. (Color online) Contour plot of the ratio 7 of elastic to 
inelastic cross sections for Li+NH collisions as a function of 
energy and magnetic field, obtained from MQDT calculations 
with perturbative long-range corrections, using a basis set 
with n max = 10 and L max = 8. The diagonal blue line shows 
the field at which the Zeeman energy in a quadrupole trap is 
6k B T. 



We have been able to reproduce the results of coupled- 
channel calculations across the entire range of energies 
and field relevant to sympathetic cooling of NH by Li, 
using our new version of MQDT combined with perturba- 
tive corrections for long-range inelasticity caused by the 
magnetic spin dipolar interaction. The MQDT results 
required coupled-channels calculations at only 5 combi- 
nations of energy and field, whereas the coupled-channels 
calculations [17] required 204 combinations. MQDT thus 
has enormous potential as an efficient computational 
method for molecular collisions. 
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